
%!TEX program = xelatex
%!TEX TS-program = xelatex
%!TEX encoding = UTF-8 Unicode

\documentclass[10pt]{article} 

\input{wang_preamble.tex}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage{titling}
\setlength{\droptitle}{-2cm}   % This is your set screw

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%文档的题目、作者与日期
%%\author{王立庆（2019级数学与应用数学1班）}
%\author{学号 \underline{\hspace{4cm}}\,\,\,\, 姓名 \underline{\hspace{4cm}}  }
%%\title{高等代数第六章：向量空间}
%\title{统计软件考试解答 }
%%\date{\vspace{-3ex}}
%\renewcommand{\today}{\number\year \,年 \number\month \,月 \number\day \,日}
%\date{2023年4月24日}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{document}

%\maketitle

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{center}

{\Large\bf \H 上海立信会计金融学院期终考试卷 } \hspace{0.3cm} {\Large \underline{ A }卷 解答}

\vspace{0.3cm}

{\large \bf \H 2023 $\sim$ 2024 学年 第 二 学期 }

\vspace{0.3cm}

{\large \bf \H \underline{ \emph{2021-2022级数学与应用数学专业} } 《\underline{ \emph{应用数学前沿专题} }》 课程代码：\underline{ 16330012B }  }

\end{center}

\vspace{0.3cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{enumerate}

%1:素数的筛法。 
%2:Mandelbrot集合的绘图。 
%3:Julia集合的绘图。 
%4:符号计算隐函数绘图。 
%5:Van de Pol方程求解。 
%6:Sturm-Liouville方程求解。 
%7:Bratu方程求解。 
%8:延迟微分方程求解。 
%9:使用离散余弦变换压缩图像。 
%10:使用神经网络模型识别手写数字。 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %1
\begin{enumerate}
\item  Eratosthenes筛法是如何找出所有素数的？
\item  使用下述程序计算1-1000之间有多少个素数。
\item  这段程序的while语句做了什么？
\end{enumerate}

\begin{python}
def sieve_v1(n):
    primes = list(range(2,n+1))
    for p in primes:
        if p*p>n:
            break
        product = 2*p
        while product <=n:
            if product in primes:
                primes.remove(product)
            product+=p
    return len(primes), primes
\end{python}

\vspace{0.2cm}


{\color{red}解答：
\begin{enumerate}
\item  Eratosthenes筛法把自然数从2开始，从小到大排成一列，然后保留2，划去后面的2的倍数，然后保留2之后的第一个剩下的数，也就是3，然后划去后面的3的倍数，然后保留5，然后划去后面的5的倍数，如此重复。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\item  在1-1000之间有168个素数，分别为 $2,3,5,7,11,\cdots, 983,991,997$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\item  这段程序的while语句把素数 $p$ 的2倍、3倍、等等数字从这个列表中划去。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\end{enumerate}

}

\vspace{0.2cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %2
称复数 $c$ 是 Mandelbrot 集合中的点，如果数列 $$z_0=c, z_{n+1}=z_n^2+c, n=0,1,2,\cdots$$ 是有界的。
设 $c$ 不是 Mandelbrot 集合中的点。设 $N$ 是使得 $|z_N(c)|>2$ 的最小的正整数。
称这个 $N$ 为 $c$ 的逃逸参数。
\begin{enumerate}
\item  参考下述程序，编程计算 $c=-1.17+0.2j$ 的逃逸参数。
\item  查阅资料或目测，Mandelbrot 集合的面积的大概是多少？
\item  证明：若 $|z_n|\ge 2$, 则 $|z_{n+1}|\ge 2$. 
\end{enumerate}

%编写一个循环语句，每次判断这个数列的下一项的模长是否大于2.
\begin{python}
c=-1.17+0.2j; z=c; count=0
while (abs(z)<=2) and (count<1000):
    z=z**2+c
    count=count+1
print(count)
\end{python}

%另一种计算，直接写出这个数列。
%\begin{python}
%#c=-1.16+0.2j
%c=-1.17+0.2j; z=c; x=[]; y=[]
%for k in range(100):
%    x.append(z.real); y.append(z.imag)
%    if abs(z)<2:
%        print('%d : %5.2f + %5.2f j , norm: %5.2f.' %(k,z.real,z.imag,abs(z)) )
%    else:
%        break
%    z=z**2+c
%print(z)
%import matplotlib.pyplot as plt
%plt.plot(x,y,'.-')
%\end{python}

\vspace{0.2cm}


{\color{red}解答：

\begin{enumerate}
\item  程序计算可知，当 $N=27$ 时， $|z_N(c)|=1.30$, 而当 $N=28$ 时， $|z_N(c)|=2.78$. 
因此这个 $c=-1.17+0.2j$ 的逃逸参数是 $N=28$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\item  Mandelbrot 集合的面积大概比1.5多一些。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\item   设 $|z_n|\ge 2$. 根据复数模长的三角不等式，可得
$$
|z_{n+1}| = |z_n^2+c| \ge |z_n|^2 - |c| \ge 4 - 2 = 2, 
$$
其中直接计算验证可得 $|c| = \sqrt{(-1.17)^2+(0.2)^2} < 2$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

%\url{https://people.math.harvard.edu/~knill/teaching/math21a2019/exhibits/mandelbrot/mandelbrot.pdf} 

\end{enumerate}

\begin{figure}[ht]
\centering
\includegraphics[height=4cm, width=6cm]{final-a-p2-mandelbrot.png}
\caption{第2题：判断 Mandelbrot 集合附近的一个点生成的数列是否有界}
\end{figure}



}

\vspace{0.2cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %3

考虑范德波尔方程 $y''-(1-y^2)y'+y=0$, 设初始条件为 $y(0)=2, y'(0)=3$. 
\begin{enumerate}
\item  将该方程写成标准的一阶微分方程组。
\item  使用 scipy.integrate 模块的 odeint 函数计算数值解，求出 $y(10)$ 和 $y'(10)$ 的值。

\begin{python}
import numpy as np
from scipy.integrate import odeint

def rhs(z,t):
    return [ z[1], (1-z[0]**2)*z[1]-z[0] ]
t=np.linspace(0,10,101)
y0=np.array([2,3])
y=odeint(rhs,y0,t)
import matplotlib.pyplot as plt
plt.plot(y[:,0],y[:,1],'.-')
\end{python}

\end{enumerate}

\vspace{0.2cm}



{\color{red}解答：
\begin{enumerate}
\item  将未知函数 $y(t)$ 与其一阶导函数 $v(t):=y'(t)$ 写成一个向量 $\vec{z}(t)=(y(t),v(t))^T$, 
%begin{eqnarray*}
%\vec{z}(t) = \begin{bmatrix} y(t) \\ v(t) \end{bmatrix}, 
%\end{eqnarray*}
则范德波尔方程的标准形式是 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\begin{eqnarray*}
\frac{d}{dt}{\vec{z}} (t) = \vec{f}(\vec{z}), t\ge 0, \,\,\, \vec{z}(0)=\vec{z}_0, 
\end{eqnarray*}
其中方程的右端函数和初值条件分别为
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\begin{eqnarray*}
\vec{f}(\vec{z}) = \begin{bmatrix} v \\  (1-y^2)v-y \end{bmatrix}, \hspace{0.3cm}
\vec{z}_0 = \begin{bmatrix} y_0 \\ v_0  \end{bmatrix} = \begin{bmatrix} 2 \\ 3  \end{bmatrix}. 
\end{eqnarray*}

\item  使用程序计算可得 $y(10)=0.1605$ 和 $y'(10)=-2.0074$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\begin{figure}[ht]
\centering
\includegraphics[height=4cm, width=6cm]{final-a-p3-van-de-pol.png}
\caption{第3题：范德波尔方程}
\end{figure}

\end{enumerate}

}

\vspace{0.2cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %4
\begin{enumerate}
\item  麦克-格拉斯方程是生物数学里的一个延迟微分方程模型，
\begin{equation*}
\left\{
\begin{aligned}
\frac{dx(t)}{dt} &= f(t,x) = a\frac{x(t-\tau)}{1+x(t-\tau)^m} - bx(t),\,\, t>0, \\ 
x(t) &= x_0(t), \,\, t\in [-\tau,0]. 
\end{aligned}
\right. 
\end{equation*}
其中 $a,b,\tau,m$ 是非负常数。
请介绍有关延迟微分方程的数值解的软件与参考文献。

\item  对 $-\tau\le t\le 0$, 设有初始值 $x(t)=x_0$. 使用一阶欧拉方法 $x_{i+1}=x_i+hf(t_i,x_i)$, 编程实现上述方程的数值解。设 $a=2,b=1,\tau=5,m=20,h=0.1,x_0=0.5,-5\le t\le 100$. 则 $x(100)$ 等于多少？

\begin{python}
import numpy as np
import matplotlib.pyplot as plt
a=2; b=1; tau=5; m=20; x_initial=0.5; t_final=100; h=0.1; d=50 
N=(tau+t_final)*10+1; t=np.linspace(-tau,t_final,N) 
x=np.zeros_like(t); x[0:d+1]=x_initial 

s=np.zeros_like(t) 
for k in range(d,len(t)-1):
    s[k]=a*x[k-d]/(1+x[k-d]**m)-b*x[k]
    x[k+1]=x[k]+h*s[k]

fig = plt.figure()
ax = fig.add_subplot(121)
ax.plot(t,x,'b--')
bx = fig.add_subplot(122)
bx.plot(x[0:N-15],x[15:N],'--')
fig.subplots_adjust(hspace=0.5, wspace=0.5)
\end{python}

\end{enumerate}

\vspace{0.2cm}


{\color{red}解答：
\begin{enumerate}
\item  延迟微分方程的数值求解，开源软件有 python 的 pydelay, 由 V. Flunkert 与 E. Sch{\"o}ll 在2009年开发，
项目地址为 \url{https://pydelay.sourceforge.net}. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

pydelay的项目论文为 \url{https://arxiv.org/abs/0911.1633}. 
其它专著有 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
\begin{itemize}
\item  Murray, J.D., Mathematical biology I. An Introduction, Springer, 2002.
\item  Erneux, T., Applied Delay Differential Equations, Springer, 2009. 
\end{itemize}


\item  运行程序可得 $x(100)=0.04889$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\begin{figure}[ht]
\centering
\includegraphics[height=4cm, width=10cm]{final-a-p4-mackey-glass.png}
\caption{第4题：麦克格拉斯方程}
\end{figure}

\end{enumerate}

}

\vspace{0.2cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %5
记变量代换 $\left\{\begin{array}{rcl} u &=& x^2+y^2 \\ v &=& xy \end{array}\right.$ 的雅可比矩阵为 $J(x,y)$. 
\begin{enumerate}
\item  求雅可比行列式 $|J(x,y)|$ 在点 $(x,y)=(1,2)$ 的值。
\item  雅克比矩阵的几何意义是什么？
\end{enumerate}

\vspace{0.2cm}


{\color{red}解答：

\begin{enumerate}
\item  按照雅克比矩阵的定义，计算可得 $J(x,y)=\begin{bmatrix} u_x & u_y \\ v_x & v_y \end{bmatrix} = \begin{bmatrix} 2x & 2y \\ y & x \end{bmatrix}$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}}) 

计算行列式可得 $|J(x,y)|=2x^2-2y^2$.代入点可得 $|J(1,2)|=-6$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\item  多元向量值函数在局部范围内可以近似为一个线性映射，雅克比矩阵是这个线性映射的系数矩阵。
在一元函数的特殊情形，雅克比矩阵就是导数，也就是切线的斜率。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\end{enumerate}

}

\vspace{0.2cm}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %6
考虑伯格斯方程的初边值问题 
\begin{equation*}
\left\{ 
\begin{aligned}
& u_t = -uu_x +\mu u_{xx}, \,\, 0\le t\le T, a\le x\le b, \\  
& u(0,x) = f(x), \\ 
& u(t,a) = g_1(t), u(t,b) = g_2(t).
\end{aligned}
\right. 
\end{equation*}

\begin{enumerate}
\item  有限差分法的求解思路是什么？
\item  谱方法的求解思路是什么？
\item  物理信息神经网络(PINN)的求解思路是什么？
\end{enumerate}

\vspace{0.2cm}


{\color{red}解答：

\begin{enumerate}
\item  有限差分法将时间和空间划分为一系列网格点，用差分来近似微分，将微分方程化为代数方程求解。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\item  谱方法将解函数写成一些正交基函数的组合，组合系数待定。通过傅立叶变换将微分方程化为代数方程，求出待定系数。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\item  神经网络的输入为时间和空间坐标，输出为目标物理量，如伯格斯方程的速度场。用初边值条件和微分方程来定义损失函数，使用反向传播和梯度下降等优化方法，最小化整个损失函数，则输出结果将尽可能地符合微分方程和初边值条件。
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\end{enumerate}

}

\vspace{0.2cm}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage 
\item %7

长度为 $N$ 的一维数组 $(f(0),f(1),\cdots,f(N-1))$ 的离散余弦变换(DCT)是指另一个长度为 $N$ 的一维数组
 $(g(0),g(1),\cdots,g(N-1))$, 其中每一个分量为
\begin{eqnarray*}
g(t) = a(t) \sum\limits_{x=0}^{N-1}f(x)\cos \left[ \frac{\pi(2x+1)t}{2N} \right], \,\,\,\, t=0,1,2,\cdots,N-1, 
\end{eqnarray*}
其中 $a(0)=\sqrt{1/N}$, $a(1)=a(2)=\cdots=a(N-1)=\sqrt{2/N}$.
\begin{enumerate}
\item  写出 $N=4$ 时的变换矩阵。
\item  写出傅立叶变换的定义。
\end{enumerate}

\vspace{0.2cm}


{\color{red}解答：

\begin{enumerate}
\item  从 $f$ 到 $g$ 的变换系数组成一个正交矩阵，
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})
{\small  
\begin{eqnarray*}
\begin{bmatrix} g(0) \\ g(1) \\ g(2) \\ g(3) \end{bmatrix}
=
\begin{bmatrix}  \sqrt{1/4}& \sqrt{1/4} & \sqrt{1/4} & \sqrt{1/4} \\  
\sqrt{2/4}\cos(\pi/8) & \sqrt{2/4}\cos(3\pi/8) & \sqrt{2/4}\cos(5\pi/8) & \sqrt{2/4}\cos(7\pi/8)  \\  
\sqrt{2/4}\cos(2\pi/8) & \sqrt{2/4}\cos(6\pi/8) & \sqrt{2/4}\cos(10\pi/8) & \sqrt{2/4}\cos(14\pi/8)  \\  
\sqrt{2/4}\cos(3\pi/8) & \sqrt{2/4}\cos(9\pi/8) & \sqrt{2/4}\cos(15\pi/8) & \sqrt{2/4}\cos(21\pi/8)  \\  
\end{bmatrix}
\begin{bmatrix} f(0) \\ f(1) \\ f(2) \\ f(3) \end{bmatrix}.
\end{eqnarray*}
}

\item  设 $f:\mathbb{R}\to\mathbb{R}$ 是一个定义在实数轴上的实值函数，设 $f$ 是绝对可积的。则 $f$ 的傅立叶变换是指函数 $g: \mathbb{R}\to\mathbb{C}$, 由下式定义
\begin{eqnarray*}
g(t) = \frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty} f(x)e^{-itx}dx.
\end{eqnarray*}
函数 $f(x)$ 的傅立叶变换 $g(t)$ 经常记为 $\mathcal{F}[f](t)$ 或者 $\hat{f}(t)$. 
\dotfill (\underline{\hspace{0.2cm} 5分 \hspace{0.2cm}})

\end{enumerate}

}

\vspace{0.2cm}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\end{enumerate}

\end{document}





